% DETERMINE WHICH TRIALS TO OMIT FROM GRAND AVERAGES
% for slow speed data (18 deg/s)

omit = [];
p = 1;
for t = 1:numel(files)
    if trialsummary(1,t).RMSE > outlier  % adjust RMSE cut off if nec
        omit(p) = t;
        p = p + 1;
    end
end

% CALCULATE AVERAGE TRACES & RMSE'S OVER ALL VALID TRIALS (3 CYCLES)

sn = 0;
ss = 0;
Empty = zeros(18000,1);

SN.X = Empty;
SS.X = Empty;

SN.allerrsq = [];
SS.allerrsq = [];

SN.allerrsqocc = [];
SS.allerrsqocc = [];

SN.allerrsqfix = [];
SS.allerrsqfix = [];

for t = 1:numel(files)
    if ismember(t,omit) == 0  % is trial valid?
        
        if trialsummary(1,t).stimulus.smooth == 0
            SN.X(:,sn+1) = trialsummary(1,t).data.X;
            SN.velnosac(:,sn+1) = trialsummary(1,t).data.vel;
            SN.allerrsq = [SN.allerrsq;...
                trialsummary(1,t).errsq];
            SN.allerrsqocc = [SN.allerrsqocc;...
                trialsummary(1,t).errsqocc];
            SN.allerrsqfix = [SN.allerrsqfix;...
                trialsummary(1,t).errsqfix];
            sn = sn + 1;
        else
            SS.X(:,ss+1) = trialsummary(1,t).data.X;
            SS.velnosac(:,ss+1) = trialsummary(1,t).data.vel;
            SS.allerrsq = [SS.allerrsq;...
                trialsummary(1,t).errsq];
            SS.allerrsqocc = [SS.allerrsqocc;...
                trialsummary(1,t).errsqocc];
            SS.allerrsqfix = [SS.allerrsqfix;...
                trialsummary(1,t).errsqfix];
            ss = ss + 1;
        end
    end
end

SN.mean = nanmean(SN.X,2);
SN.mean = SN.mean(1:msPerCycle*3);
SS.mean = nanmean(SS.X,2);
SS.mean = SS.mean(1:msPerCycle*3);

SN.total = sn;
SS.total = ss;

SN.time = 1:msPerCycle*3;
SS.time = 1:msPerCycle*3;

SN.RMSE = sqrt(nanmean(SN.allerrsq));
SS.RMSE = sqrt(nanmean(SS.allerrsq));
Stats.RMSE(:,1) = [SN.RMSE, SS.RMSE];

SN.RMSEocc = sqrt(nanmean(SN.allerrsqocc));
SS.RMSEocc = sqrt(nanmean(SS.allerrsqocc));
Stats.RMSE(:,2) = [SN.RMSEocc, SS.RMSEocc];

SN.RMSEfix = sqrt(nanmean(SN.allerrsqfix));
SS.RMSEfix = sqrt(nanmean(SS.allerrsqfix));
Stats.RMSE(:,3) = [SN.RMSEfix, SS.RMSEfix];

% CALCULATE AVERAGE OVER ALL CYCLES

SN.one(:,1) = SN.mean(1:msPerCycle,1); 
SN.one(:,2) = SN.mean(msPerCycle+1:msPerCycle*2,1);
SN.one(:,3) = SN.mean(msPerCycle*2+1:msPerCycle*3,1);
SN.one = nanmean(SN.one,2);

SS.one(:,1) = SS.mean(1:msPerCycle,1); 
SS.one(:,2) = SS.mean(msPerCycle+1:msPerCycle*2,1);
SS.one(:,3) = SS.mean(msPerCycle*2+1:msPerCycle*3,1);
SS.one = nanmean(SS.one,2);
      
% CALCULATE MEAN VELOCITIES EXCL SACCADES & SMOOTH

SN.meanvelnosac = nanmean(SN.velnosac,2);
SN.meanvelnosac = smooth(SN.meanvelnosac(1:msPerCycle*3),smve);
SS.meanvelnosac = nanmean(SS.velnosac,2);
SS.meanvelnosac = smooth(SS.meanvelnosac(1:msPerCycle*3),smve);

SN.onevelnosac = (SN.meanvelnosac(1:msPerCycle,1) + ...
                  SN.meanvelnosac(msPerCycle+1:msPerCycle*2,1) + ...
                  SN.meanvelnosac(msPerCycle*2+1:msPerCycle*3,1))./3;
              
SS.onevelnosac = (SS.meanvelnosac(1:msPerCycle,1) + ...
                  SS.meanvelnosac(msPerCycle+1:msPerCycle*2,1) + ...
                  SS.meanvelnosac(msPerCycle*2+1:msPerCycle*3,1))./3;

% CALCULATE MEAN VELOCITIES INCL SACCADES (smoothing sm = 100)

SN.meanvel = velocity(SN.mean,smvi,scnw,scnwcm,view_dist);
SS.meanvel = velocity(SS.mean,smvi,scnw,scnwcm,view_dist);

SN.meanvel = SN.meanvel(1:msPerCycle*3);
SS.meanvel = SS.meanvel(1:msPerCycle*3);

SN.onevel = velocity(SN.one,smvi,scnw,scnwcm,view_dist);
SS.onevel = velocity(SS.one,smvi,scnw,scnwcm,view_dist);

SN.onevel = SN.onevel(1:msPerCycle,1);
SS.onevel = SS.onevel(1:msPerCycle,1);

% CALCULATE LAGS

msPerhCycle = floor(msPerCycle/2);

SN.lag(1:msPerhCycle,1) = target18.X(1:msPerhCycle,1) - SN.one(1:msPerhCycle,1);
SN.lag(msPerhCycle+1:msPerCycle,1) = SN.one(msPerhCycle+1:msPerCycle,1) - ...
    target18.X(msPerhCycle+1:msPerCycle,1);
SS.lag(1:msPerhCycle,1) = target18.X(1:msPerhCycle,1) - SS.one(1:msPerhCycle,1);
SS.lag(msPerhCycle+1:msPerCycle,1) = SS.one(msPerhCycle+1:msPerCycle,1) - ...
    target18.X(msPerhCycle+1:msPerCycle,1);